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Abstract 

A general method for calculating next-to-leading order cross sections in perturbative 
QCD is presented. The algorithm is worked out for calculating iV-jet cross sections 
in hadron-hadron collisions. The generalization of the scheme to performing caclula- 
tions for other QCD process, such as electron-positron annihilation or in deep inelastic 
scattering is also straightforward. As an illustration several three-jet cross section 
distributions in electron-positron annihilation, calculated using the algorithm, are pre- 
sented. 
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1 Introduction 



It is well known that owing to large scale dependence any result of a leading order calculation 
in perturbative QCD can be regarded only as a one parameter fit to the data, but not a real 
theoretical prediction. The scale dependence is expected to decrease substantially with the 
inclusion of next-to-leading order corrections, which in a sense fixes the scale. Consequently, 
the calculation of cross sections at next-to-leading order in perturbative QCD is highly desir- 
able. During the last sixteen years a good number of such calculations have been performed. 
A common feature of these calculations is that they use the specific, usually simple kine- 
matics of the problem in order to achieve the analytic cancelation of infrared singularities. 
This specialization brings an element of art into next-to-leading order calculations, which 
for more complex problems like the calculation of four-jet production in e + e~ annihilation 
or three-jet production in hadron-hadron collisions becomes an obstacle to overcome. One 
would rather like to apply a sort of "standard technique" in order to discriminate the prob- 
lems of book-keeping from those of the theory. The first steps into the direction of drawing 
a standard picture were done by Ellis, Kunszt and Soper in refs. [|T], and by Giele, Glover 
and Kosower in refs. || |J. The common characteristics of these works is the recognition 
that the factorization properties of the QCD amplitudes or squared matrix elements |J as 
well as that of the phase space in the limits when one particle becomes soft, or an external 
pair becomes collinear can be utilized for devising a universal scheme for the calculation of 
any infrared safe physical cross section at next-to-leading order in perturbative QCD.[] 

The algorithm developed in refs. fl], @] - - called the subtraction method — was applied 
for calculating inclusive one-jet |7fl as well as for two-jet [|J production at next-to-leading 
order. The same algorithm however, cannot be directly applied to more complex cases 
- like the ones mentioned above — , because the evaluation of certain integrals used the 
specific 2^2 kinematics of the problem considered and also because the algorithm relies 
on single singular decomposition of the squared matrix element. The latter is not a problem 
in principle, however, the last decade has proved that one has to use helicity amplitudes 
both at tree (see eg. M) as well as at loop level [ffDl O, O in order to obtain higher order 
results and it is rather cumbersome to square these amplitudes analytically and perform the 
single singular factoring of the squared matrix elements. The subtraction method has been 
generalized to the calculation of three-jet cross sections in hadron collisions in a recent paper 



by Frixione, Kunszt and Signer ||13|| . In this paper however, the physical quantity — the 
"measurement function" — was used for coping with the problem of single singular factoring, 
which makes the generalization to different types of cross sections than the one discussed in 
the paper non-trivial. 

The algorithm of refs. || || - called the slicing method — avoid the above obstacles 
offering a general scheme, but at the price of introducing an unphysical parameter s m i n and 
calculating the result to 0(s m ; n ) accuracy. In principle s m i n can be chosen infinitesimal, 
thus an exact result can be recovered. However, in practice the choice of a very small 
s m i n adversely affects the numerical convergence of the Monte Carlo integrals and one has 
to carry out a balancing procedure between the error of the Monte Carlo integration and 
the one introduced by the choice of finite s m i n in order to minimize the theoretical error. 



^^More recently an elegant scheme has been outlined by Catani and Seymour in ref. 
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This balancing procedure can be inconvenient in those cases when the matrix elements are 
complicated and their numerical evaluation is time consuming. 

The aim of the present paper is to provide a simple generalization of the subtraction 
method that can be used for calculating any infrared safe physical quantity at next-to- 
leading order in perturbative QCD if the required tree and one-loop level helicity amplitudes 
are known. In order to minimize the theoretical error of the Monte Carlo integrals, we apply 
important sampling. There are two ways of achieving efficient important sampling. One is 
when the integrand is decomposed into single singular factors and the important sampling 
is performed in the variable controlling the singularity. As stated above, however, single 
singular factoring is better avoided. The second possibility, that we apply in this paper, is a 
decomposition of the phase space into regions, where the integrand can become singular due 
to the vanishing of only one Lorentz invariant of the external momenta. This decomposition 
can be done quite generally, without any reference to the squared matrix element, or to 
the physical quantity being calculated. We describe the algorithm in detail for the case of 
hadron collisions, which is the most general case one can encounter. Algorithms for other 
processes can be obtained by leaving out certain terms as it will be explained later. 

In section 2, we discuss how infrared safe cross sections can be calculated in perturbative 
QCD. In the following sections we describe the cancelation scheme in detail. The scheme 
is based upon the soft and collinear factorization properties of the squared matrix elements 
of QCD and that of the phase space. The singularity structure of the one-loop amplitudes 
for QCD processes involving arbitrary number of external partons has been discussed in ref. 
l3| . Using those results it is not difficult to find the universal structure of the singularities in 



the next-to-leading order matrix element of the virtual corrections. That universal structure 



has already been given in refs. |2|, ITo] . For the sake of completeness as well as for setting 
some of the notation, we recall the necessary formulas in section 3. 

Section 4 contains the essence of our algorithm. Here we discuss the decomposition of 
the phase space, the singularity structure of the real corrections. We describe how the phase 
space is generated to achieve the necessary important sampling. We define local soft and 
collinear subtraction terms that make the integral of the real corrections over the N + 1 
particle phase space finite. The explicit expression for this finite integral is also given. We 
integrate out the variables of the soft or collinear particle analytically in sections 5 and 6. We 
show that the remaining expression has the form of the 2 — > N integrals (like the Born and 
virtual corrections), so they can be combined and the analytic cancelation of the infrared 
divergences is demonstrated. The remaining finite 2 — > iV integral is explicitly given. Section 
7 contains some sample results for three-jet cross section calculation in e + e~ annihilation. 
We conclude in section 8. The appendix is a collection of the analytic integrals that were 
used in the main text for the demonstration of the cancelation of the infrared divergences. 



2 Infrared safe cross sections at next-to-leading order 

At order ay +l \ one calculates cross sections, with infrared divergences controlled using 
dimensional regularization in d = 4 — 2e dimensions, for the two incoming hadrons to collide 
and produce either N or N + 1 final state partons. According to the factorization theorem, 
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the next-to-leading order infrared safe physical cross section is a sum of two integrals, 



where 



a = J[2 -> N] + J[2 -> N + 1], (2.1) 



J[2->W] = X! J dx A f A (a A ,x A ) J dx B f B (a B ,x B ) (2.2) 



and 



/ dr^W, . . . ,p&)(|A<(2 - iV)| 2 )^(rf, . . . 
J[2 -> iV + 1] = £ / dx A / A (a A , x A ) / dx B f B (a, B , x B ) (2.3) 

a A ,a s ,ai,...,a JV+ i 

i y dr<" +1 >(p?, . . .,p% +l )(\M(2 1)\ 2 )S N+ M, . . . ,p% +1 ). 
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In these equations s = x A x B s. In the phase space measures, 

dK->«, . . . , K) = 1 ft (gp^) <*0 V* (rf + Pb - Erf) (2-4) 

(n = iV, A" + 1) we included the identical particle factor l/n\ which is present if we treat 
all final state partons identical and sum over all possible parton flavors Oj = u,u,d,d, . . . , g. 
The parton distribution functions for incoming partons A and B defined in the MS renor- 
malization scheme are denoted by f A (a A ,x A ) and f B (a B ,x B ). 

In order to factor the dependence of the cross section on the physics of low transverse 
momenta out of the partonic cross section and into these MS parton distributions, in the 
2 — > iV cross section one uses the modified parton distribution f(a,x) that satisfies 

f(a,x) = Y: f—f(b,x/z) (2.5) 



b 

X 



<W ~z) + J*^' ^Pa /b (z) + O^) 



with P a /b(z) being the full Altarelli-Parisi kernel, for the b — * a splitting: 

P a/b (z) = P a/b (z) - 5 ab 2 -^ + S ab ^ C{a ] + 5 abl (a)5(l - z), (2.6) 

l — z [i — z) + 

where, for instance, in the case of g — > gg splitting 

P g/g (z) = 2C(g) + _L_ + *(1 - *)) ■ ( 2 - 7 ) 

C{g) = N c is the color charge of the gluon, while C(q) = V/(2N C ) is that of the quark 
(V = (iV c 2 — 1)), while the 7 constants represent the contribution of the virtual graphs to 
the Altarelli-Parisi kernel, 

llN c -2N f W , . 

7( 9 ) = ^j-^, 7(9) = j^- (2.8) 
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The notation 1/(1 — z) + is the usual "+" prescription, 



i f(z) _ M{z)-f{l) 

o (l-z) + Jo l-z 1 ' ' 



This factorization recipe is discussed in ref. fli~6 |. The (|-M(2 — > n)\ 2 ) functions are the 2 — > n 
squared matrix elements averaged over initial state and summed over final state spins and 
colors: 

(\M(a + b^n)\ 2 ) = ^— £ \M(a + b^n)\\ (2.10) 

\ color/ 

In the conventional MS scheme, we need (|7V1| 2 ) in d = 4 — 2e dimensions. However, it was 
shown in ref. [[DJ that simple rules exist which tell us how to obtain the finite 2 — > N hard 
scattering cross section of the conventional MS scheme at next-to-leading order using the 
expressions for (|.M(2 — > n)\ 2 ) obtained in the dimensional reduction scheme. Therefore, we 
use u(g) = 2V and w(q) = 2N C , which are valid in d — 4 dimensions and the four- dimensional 
expressions for the squared matrix elements. Finally, the functions S n define the physical 
quantity to be calculated. 



In equation (|2.1|) both terms are singular when the regularization is removed, e — > 0. 
When e ^ the singularities are represented as 1/e 2 and 1/e poles. These poles cancel 
between the X[2 — > N] and X[2 — > N + 1] terms, provided the physical measurement, rep- 
resented by the functions «S n , is infrared safe. This means that the emission of a soft or a 
collinear parton must not influence the result of the measurement. Therefore, the measure- 
ment functions S n must possess the following properties: 

lim S N+1 (p^,...,p^ +1 ) =<Sjv(Pi, >£, ...,p% +1 ), ie[l,JV + l]; (2.11) 



^■p^ 1 — *■(! — z)p 



i,je[i,N + l], ze[o,l], 

(2.12) 

lim S N+1 (p1, . . . ,p% +1 ) = S N {pt, .,p% +1 ), 
m = A,B, je[l,N + l], z€[0,l]. 



3 The 2 -> N integral 

In this section, our aim is to write the X[2 — > N] integral in such a form that will make 
the cancelation of divergent pieces against corresponding divergent terms in the X[2 — > 
N + 1] integral as simple as possible. The discussion is a generalization of the corresponding 
discussion in ref. @ given for the X[2 — ► 2] integral to the 2 —>■ N case. There are some 
differences however. Firstly, we do not specify the integration variables, but leave it for the 
reader to use a preferred choice. Secondly, in order we could use the results for the five-parton 
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one-loop QCD helicity amplitudes [K], [LI], [12j for a next-to-leading order calculation of three- 
jet production in hadron collisions, we perform the analysis using dimensional reduction 
scheme and add the necessary transition terms at the end to obtain the correct formula in 



conventional MS scheme as they are given in ref. [15| . 

In order to simplify the book-keeping of the various factors, we introduce the integration 
measure 

1 1 n N N 
V N {e) = YsMj2^Z dXAdXB I! [(2vr/i) 2£ d 4 - 2 ^25( ft 2 )] , (3.1) 

where the subscript on V reminds us that this measure is related to the phase space inte- 
gration measure of iV particles according to the relation 

{ ^f^dx A dx B dT^(pl ...,p%)= V N (e)(2n^6 d U+Pb ~ fx) • ( 3 - 2 ) 

It is convenient to write the perturbative expansion of the squared matrix element 
summed over final spins and colors and averaged over initial spins and colors and with 
ultraviolet renormalization in the MS renormalization scheme included in terms of functions 
and *g? +2) (p), 

(\M(2^N)\ 2 ) ^ (3.3) 

where 

and the subscript DR refers to expressions obtained using dimensional reduction. We use 
the notation a = (a a, clb, ■ ■ ■ , %) for the collection of external parton flavors and p = 
{PaiPbiPit ■ ■ iPn) f° r the collection of external four-momenta. The variable Qes is an 
arbitrary parameter of mass dimensions introduced to facilitate writing the result [pL7| . The 
dependence of the function ty^R + ^ on Qes is such that the squared matrix element does not 
actually depend on Qes- The first term in the curly braces is the Born 2 — > N matrix element 
squared, without the g 2N coupling factor, while the second term is the next-to-leading order 
contribution, which in the dimensional reduction scheme can be expressed in terms of the 
helicity amplitudes |15|| , 




,2 



(W^or h£r ^ +2 \a,p) = EE(^ (1) ^ (0) *+^ (1) *^ (0) ) • (3-5) 

Z7r VVES/ hel col 

On the right hand side the superscript (0) and (1) refers to the tree and one- loop helicity 
amplitudes respectively. The helicity amplitudes can be decomposed in color space in terms 
of gauge invariant color subamplitudes. For instance, in the case of pure gluon processes, 

A® (9i, ■■■,9n)=9 N (£) E C^»*»aW(l, . . . , N), (3.6) 
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where C 9u '" ,9N is a general invariant matrix in color space with upper indices in the adjoint 
representation and summation on n runs over a linearly independent set of such matrices. 
Explicit examples of this decomposition can be found in reference [0. The one- loop color 
subamplitudes and thus the helicity amplitudes can naturally be decomposed into singular 
terms containing at most double poles in e and into terms that are finite when e — > 0, 



(3.7) 



Looking at the explicit form of the singular terms of one-loop five-parton color subamplitudes 
n], we see that the imaginary parts of the factors — l/e 2 (— Sij/Q^ s )~ £ do not contribute to 
the function ^dr. In such cases the ^dr functions have the following structure: 



(2N+2) 
DR 1 



a,p) 



1 



1J1; K a,p) \ -— E C ( l 

E n=A,BA,...,N 



E iM 



(3.8) 



n=A,BA,...,N 



1 



TV 



^ c (m,n=A,B,l 1 ...,N\ 



y^mn) ^ van, DR 



TV 



E tK) -7 E 

n=A,B,l,...,N (m,n=A,B,l,...,N\ 



+2 



(4vra s ) — c r 

Z7T 



2 

ES, 



VE(^ (o) *+4ri (o) )+°( £ )- 



hel col 



Here the \E f ^ ;c - ) functions are the color correlated Born squared matrix elements defined 



in ref. [HI. The factor 2 



(47r« s )^(« s /27r)c r (/i 2 /g| s ) £ 



S 2Af (5/47T) 



Ofe) cancels 



against the coupling factors in A^gA/ - 1 . The functions £(x) and ^(^) are defined as 

x 



£(x) = In 



Q 2 



ES 



£ 2 (x) =t(x)-7r'0(x). 



(3.9) 



Substituting the integration measure of eq. fl3.2|), the perturbative expression for the 
squared matrix element, eq. ( [3.3|) and the expression for the modified effective parton dis- 
tribution functions as defined in eq. ( |2.5p , we can write the 2 — > N cross section as 

J[2 -> iV] DR = (3.10) 
E / V N {e)S N {rt, . . .,p» N )(2n^5 d (p A +p B - Erf) 



a A ,a B 



x jL(a4, a B ,x A ,x B 
w(a' A ) f 1 dz 



(27V)/- — ** . ft s 

dr fl>PJ + 7-cr 

Z7T 



,Qes, 
(4tt) £ 



(27V+2) 
DR 1 



a,p) 



o;(a^) Z" 1 dz 



—zLi a A , a B , x A , — 



e) 2tt 
(4tt) £ a 



DR 1 



z J elYl-e) 2 



(27V) 
DR 1 



a,p) 



a,p) \. 
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The function L used here and elsewhere describes the parton luminosity: 

Tf \ f(a A ,x A ) f(a B ,x B ) , Q11 ^ 

L{a A ,a B ,XA,XB) = —, — ; — r • (3-H) 

bj{a A )XA uj{a B )x B 

According to the factorization and Kinoshita-Lee-Nauenberg theorems, and we shall see it 
at next-to- leading order explicitly in the following sections, the pole terms in eq. ( |3.10 ) cancel 
against poles emerging in the phase space integral of the brehmsstrahlung contributions. 
Therefore, it is only the Born function \I/( 2Ar ) and the finite part of the \I/( 2Ar+2 ) function - 
the last two lines of eq. ( |3.8| ) — that is really integrated in eq. ( |3.10| ). However, we need 
the corresponding finite expressions valid in conventional dimensional regularization. In ref. 
[[TJJ it was shown that simple terms are to be added in order to obtain the correct formula 
we need for a next-to-leading order calculation in the conventional MS scheme. Thus the 
function resulting from the loop corrections that we need for a next-to-leading order Monte 
Carlo program is the non-singular function 



(2iV+2) 
NS 1 



1 

4 



E 



E [% 2 hK) - 7M] + 

n=A,B,l,...,N 
,(27V;c) 



(3.12) 



'm,n = A,B,l,...,N\ 



+ 



,2AT 



47T 



-2(Smn)*^DH.(£>£) 



hel col 



where the transition terms j(a n ) are given by 



1(9) = \C{ 9 ), 7(g) = \c{q). 



(3.13) 



4 The 2 -> N + 1 integral 

In this section, we separate the X[2 — * + 1] integral into terms containing 1/e 2 and l/e 
poles, which cancel against the corresponding poles of the X[2 — > N] integral, and terms that 
are finite when e — > and, therefore, can be integrated numerically. The l/e p singularities 
arise from integrating the square of the matrix element over the (N + l)-particle phase 
space when a gluon becomes soft, or two partons become collinear. Firstly, we organize the 
integration domain so as to reduce the complexity of the problem. 



4.1 The domain of integration 

We must integrate over the momenta of the A^+l final state particles treating them identical. 
If we do not fix a definite label to each particle, then we integrate over each event topology 
(N + 1)! times. We can however, simplify the calculation by 

1. first splitting the phase space in the parton-parton cm. system into two parts: in the 
first one, the smallest angular distance = Sijf(EiEj) is between final state particles, 
while in the second region the smallest angular distance is between an initial state 
particle and a final state particle; 
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2. secondly cutting into the first region by fixing that label of the smallest Lorentz in- 
variant of final state particle pairs to which the smaller energy in the parton-parton 
cm. system belongs to be j = N + 1, and cutting into the second region by fixing the 
final state label of the smallest Lorentz invariant of pairs involving an initial state and 
a final state particle to be j — N + 1. 

With this distinction of parton (N + 1), we have to integrate over each event topology only 
N\ times and there is a corresponding symmetry factor 1/N\ associated with the integration: 

i N+l / 2eJd-l n . \ ( N+l \ 

x[0(ril +1;i) > r%£M)e{8$jP > s™ +1 )Q(E h > E N+1 ) 

T^-A'min ^ 'min J^l^min ^ b X,N+\)h 

where 

rf min } = min(di 3 - : i, j = 1, . . . ,n, i ^ j), (4.2) 

rf min = min(rf X j : X = A, B, j = 1, . . . , n), (4.3) 

d^ +1 = mm(d J<N+1 : j — 1, . . . , N), (4.4) 

dx^+i = min(d A ,N+i, d B ,N+i), (4.5) 

with d meaning either Lorentz invariant s or angular distance r. In eq. ( |4.4| ) the index k 
denotes that j for which the minimum value is assumed. 

In the cut phase space the only singularities that can occur when a single Lorentz invariant 
vanishes are 

• parton N + 1 is soft (in both regions); 

• parton iV + 1 is collinear to a final state parton i 6 [1, N] (in the first region); 

• parton iV + 1 is collinear to an initial state parton A or B (in the second region). 

We are not interested in configurations when two Lorentz invariants involving four different 
labels vanish simultaneously, because those emerge only in (JV— l)-jet configurations that 
are not considered here. 

Next, we study the singularity structure of the squared matrix element. 



4.2 Singularity structure of the squared matrix element 

In order to find the singularity structure of X[2 — > N + 1] over the cut phase space of eq. 
( |1.1| ) explicitly, it is useful to strip off the spin and color averaging and the coupling of the 
squared matrix element. We define the function 

¥^\A, B, 1, . . . , N + 1) = ^^ (\M(2 - N + 1)| 2 ), (4.6) 



S 



where the argument denotes: 

(A,B,1,...,N + 1) = (a A ,a B ,a 1 ,...a N+1 ;p l X,p^,p^...p% +1 ). (4.7) 

The singularity structure of the function \[/( 2Ar+2 ) in four dimensions can most easily be 
obtained from the factorization properties of helicity amplitudes 0. Citing only the results, 
we find for soft gluon labeled j = N + 1 

\im¥ 2N+2 \A,B,l,...,N + l) = J2 *S? ic) (ff,P)+of 



p j (m,n=A,B,l,...N\ SmjSjn \ \/ S mj y S nj ) 

(4-8) 

which is called "soft identity" in ref. B. In order to make the cancelation of the infrared 
singularities as transparent as possible, it is useful to perform single singular decomposition 
of the eikonal factor in eq. ( |4.8| ): 

2s 2s 2s 

— ■inn ! 1 1 1 1 ^ ^°mn ^ g\ 



With this decomposition and using the symmetry of the $^ ;c ' functions in the m, n indices 
we can write eq. ([4.8|) in the form 



lim ^ N ^(A,B,1,...,N+1) = £ 5 a]g ¥™ + n 2 \a,p,p>;), (4.10) 

) 



Pj fm,n=A,B,l,...N 



where 

= - rfTs ^ c \a,p). (4.11) 

In the collinear limit of two final state partons i and j — N + 1, we introduce a pseudo 
particle P with ap flavor that splits into gluons % and j: p P =Pi + Pj- The flavor a P = ai 
if Oj = g and a P = g ii ai = q, aj = q. The momentum fraction z is defined by pf = zp P . 
Then for the collinear limit of \|>( 27V+2 ) one finds: 



lim u ¥ 2N+2 \A,B,1,-..,N + 1) (4.12) 
2 



v|/gJ +2) ( 2 , P; A B, 1, . . . , N + 1) + 0(l/V%) 



where 



4 2 J +2) (z,P;A,P,l,...,iV + l) (4.13) 
= P ai/ap (^)vI/( 2 ^(AP,..., XP,..., *...) 

+ 23?e (gp^(^)$( 2Ar )(P; A, B, X- ■■)) ■ 

In this equation \1> (27V) (A, B, . . . , %P, . . . , >(...) is the \I> (2Ar) function of 2 + iV partons 
obtained from \[/( 27V + 2 )(2 — ► iV + 1) by deleting labels i and j and adding the pseudo particle 
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label P, P ai /a P {z) is the Altarelli-Parisi splitting function for the process P —>■ ij in four 
dimensions without z — 1 regulation (eq. (|2.7|) in the case of gluon splitting). The Qp^ij(z) 
functions are calculated from the tree- level splitting amplitudes, Split^ ree (i /li , j hi ) of ref. |T9|] 
according to the formula 

4gp-«w = e E c(^ J ;p)c(^ J ;p)*spiie cc (^ % j ft ospiitr(^%j /ij r (4.i4) 

where c(i,j, P) is the color matrix of the P — > zj vertex. In the case of gluon splitting 



Qp-W*) = -2C0/)«(1 - *)^f , = *(1 - *)^f, (4.15) 



while in the case of quark splitting Q q _,ij(z) = 0, which is also understood from helicity 
conservation along fermion lines. The function $( 2Ar > does not depend on the momenta 
and pj only on their sum, p P . The $ {27V) functions are calculated from the tree-level helicity 
amplitudes as the Born function \ff( 2N \ except that the summation over the helicity of parton 
P is not carried out: 

& 2N \P;A,B,...,%.. .,*...) (4.16) 

= E E A(°\A h \B h *,...,XP + , ...,x--0 

color hj\,hg,... 

xA®{A h *,B h »,...,xp-,--; *■■■)• 

In the collinear limit of a final state parton j = iV + 1 with an initial state parton A, 
we let A split into partons P and j: = p P + p^ 1 , with momentum fraction 2 defined as 
p P = zpj±, followed by an P + B — > 1, . . . N hard-scattering process. From the crossing of A 
and P in relation ( |4.13| ), for the collinear limit of \I>( 27V+2 ) one obtains: 

Hm ¥ 2N+2 \A,B,1,...,N+1) (4.17) 

Pj-*V-—")PX 

= - A(_l)/(a^)+/(ap)^g^ + 2) (l/Z) p . ^ ^ ^ N + 1} + 0(l/^T) , 

where 

/(<?) = 0, /(g) = 1. (4.18) 

We can write the right hand side of eq. ( 4.2U| ) in a more explicit form using the crossing 
relation of the Altarelli-Parisi splitting functions, 

P b/a (z) = H-l) fia)+m ^P- a/ - b (y^ (4-19) 

where a is the antiparticle of particle a. Thus we find 

_(_l)/(^)+/(ap)^gv+2) (1/ ^ p ; A, B, 1, . . . , JV + 1) (4.20) 

ujyap ) z 



_(_l)/(aA)+/(.p) 2 Sfe (Qp^.(l/^$W(P; 5, . . . , X- 
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We close the analysis of the singularity structure of \|>( 27V+2 ) with considering the limit 
when the soft gluon is also collinear with parton m. Thus we take = (1 — z)p^ n with 
z — > 1. Using the "soft-collinear" identity of ref. M, 

£ ^n' c) ip) = 2C(a m )^ N \A, B,1,...,N), (4.21) 

(n=A,B,l,...,N\ 

V n^m } 



we obtain from eq. ([4.11 ) 



1™ . „ E ^Zn 2) (P) = 2C(a m )S ajg n * ¥ 2N \p). (4.22) 



p£-»(i-*)p£, / n=ABil> ... >W x [i - z)s mj 



In the same limit, eqs. ( 4 . 1 3|) and ([4.20 ) yield 

lim Av^+ 2 )( Z; P- J A,B,1,...,N + 1) = 2C(a m )5 aj9 -—^- ¥ 2N \p) (4.23) 

S m j (1 Z)S m j 



and 



lim — *g££ a) (l/s, P',A,B,\,...,N + 1) = 2C(a A )5 a]9 —^—¥ 2N \p), (4.24) 

s Aj ~~ z ) s Aj 



in agreement with eq. ( J4.22 ) . 

The function \I/( 2Af+2 ) does not posses any other poles when parton N + 1 is soft or 
collinear to another parton. Knowing the singularity structure of the function \I>( 2iV + 2 ) ) we 
can decompose the integral X[2 — > N + 1] into three terms: 

X[2 -> iV + 1] = X[soft] + J[coll] + X[fin]. (4.25) 

The first two of these integrals possess divergences in e, but they are sufficiently simple to 
calculate the pole structure of the Laurent expansion in e around zero analytically. The 
third one is complicated, but contains at most square-root singularities over the modified 
(N + l)-particle phase space of eq. (|4. 1|) , therefore, can be integrated numerically yielding a 
finite contribution as e — > 0. We further decompose the "soft" and "collinear" contributions 
into sums of iV + 2 terms, 



X[soft] = J [soft] m , (4.26) 

m=A,B,l,...,N 

X[coll]= £ X[coll] ro , (4.27) 



m=A,BA,...,N 



(O AT I 

where X[soft] m is associated with the integral of the soft terms ^s-mn > while X[coll] m is 

( 27V 4-2) 

associated with the integral of the collinear term ^c-mN+i- We shall call these integrals 
subtraction terms for the obvious reason that subtracting them from the X[2 — > N + 1] 
integral the finite term remains. In order to define the soft and collinear subtraction terms 
precisely, we first give a decomposition of the phase space into such regions that in any one 
of them only one pair of external particles can become collinear. 
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4.3 Decomposition of the phase space integral 

In this subsection, our goal is to write the phase space integral in those variables that allow 
the most efficient Monte Carlo integration. We write the integration measure of eq. (|4.1|) as 



where 



.mf=l rrii=A,B 



(4.28) 



/ N+l \ 

ZftTCe) =V N (e)(2rr^a^ p 26(p» Pfl ) (27r^6 d U + jf B - £ rfj (4.29) 

xe(r^ 1;,) > r^ n +1;/) )e(^ > > *Wi) 

xe(min(r i)A r + i : i = 1, . . . , N, i ^ m) > r TO)Ar +i), 

/ 7V+1 \ 

V^i(e) = V N (e) (2^) 2 M 4 - 2 ^2 5(p%) (2^)-^ U + p% - £ rfj (4.30) 

xetf n +1;/) > rJl+^jec^ > ^ +1 )e(r BA+1 > r AiV+1 ), 

with p M = P%+i an d T^n{^) is defined in equation (|3.1|) . The measure 2^+4 i s defined as 
^Jv+i whh A and -B interchanged. The advantage of this decomposition of the phase space 
should be clear: in each region there is only one pair of particles that can become collinear. 
As a result, the single singular factor decomposition of the integrand is substituted by a 
(much simpler) decomposition of the phase space. 

In order to write the integration measure Pjv+i m the required form, we utilize a vector pq 
of invariant mass Q 2 that splits into the vectors pf^ and p^ = p% + i and use the mathematical 
identity 



dg 2 

2vr 



N 



V N (e) (2^)-^ U+p^-Erf 



(4.31) 



XH ( r min > r min J W l S min > Sfc.JV+lJ"^* > A J 

xG(min(r iiiV+ i : i = 1, . . . , N, i ^ m) > r m>N+1 ). 

We can use the second 8 function for integrating over the d — 1 momenta of particle m and 
over Q 2 . We obtain 



dQ^ /i^^ 2£ d^p , , , _ _ v = ^n+^ (2 )2£ ^P 

o 2vr (27r) d - 1 2E m V P; 22 v Q Fm F / E m K ^ J E 



(4.32) 

Next we choose a coordinate system in the parton-parton cm. frame which has z axis 
showing into the direction of the three-momentum pq. We denote the polar and azimuthal 
coordinates of parton N + 1 by •& and (p, so the four- vector p^ in d dimensions is 

p M = 22(1, cos sin $ cosy?, sin $ sin .. .), < i9, cp < ir. (4.33) 
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The first four components are energy, z : x and y components of the three-momentum and 
the dots mean d — 4 unspecified components. In this coordinate system 



(27T/X) 



E 



E 



\~E~ 



2e 



dE (sin^^-^d^d 1 - 2 ^ 



(4.34) 



E £ d J B(sin^)( 1 - 2£ )d^(sm^)- 2£ d^d- 2£ fi 



with f2 being the solid angle in d — 4 dimensions. We change integration variables from 
(E, cos i?) to (z, cos u;) , where z = E m / (E m + E) , so 



— £L 



(4.35) 



and cu is the angle between the three-momenta p m and p in this coordinate system, 

. 1 — z + z cos ui 
cosv = , 



(4.36) 



with p being the ratio pq/(E m + E) that can be expressed in terms of z and coscj as 



p = \Jl — 2z(l — z)(l — coso;). 



(4.37) 



In the collinear limit of particles m and N + 1, i.e. uo — > 0, the definition of z given here 
is identical to the one given in the previous subsection and therefore, extends this variable 
to non-collinear configurations naturally. We remark that this definition of z is not boost 
invariant. 

Finally, we change variable from u to Q 2 = s m ,jv+i- The necessary relation is 



coscj = 1 — 



Q 2 



2z(l - z)E 2 Q 



(4.38) 



with Eq being dependent on Q 2 . We obtain 



-7-»fin;m/ \ 



N 



1=1 



(4.39) 



x 



1 / 2 \ 2 /r, \2e 

1 /Pq\ (2irp\ 



X 



/ Q 2 Y £ ' 



PQ 



2E 2 



2 - 



Q 2 



2z(l-z)E 2 Q/ 



[z(l - z)]- £ dz 



dQ 2 



x (sin <p)- 2£ d<pd- 2£ tt 0(4^(1 -z)E 2 Q > Q 2 ) 



xe(r; n ^ 1;l) > r^)Q(s^> > s^ +1 )Q(E k > E N+1 ) 
xe(min(r iiAr+1 : i = 1, . . . , N, i ^ m) > r miN+1 ), (p E [0,7r]. 
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In this equation the Pm = Pq — p^ momentum conservation constraint is implicitly under- 
stood. We shall use the four dimensional limit of this integration measure: 



v n+i (e = 0) 



N 

■J 1 



(4.40) 



V N (e = 0)6^ U+rf-£rf 

V i=i 

^(ff) dzd Q 2d V®(Ml-z)E 2 Q >Q 2 ) 

xe(r^ n +M) > rV£M)e{8<™ > s™ +1 )Q(E k > E N+1 ) 
xe(min(r i)A r + i : i = 1, . . . , N, i ^ m) > r m>Ar +i), f £ [~^,^]- 

We now turn to the discussion of the integration measure "P]v+i • 111 order to write it in 
the required form, we imagine a 2 — > 2 scattering, A + B — ► Q + (iV + 1), followed by the 
decay of particle Q into iV particles. We write the integration measure 2?jv+i m the form 

= fp& - Erf) , (2^) 2£ d 4 - 2 ^2 5(p%) (4.41) 

V i=i / Ip^=p1+p^-p" 

Xfc) ( r min > r min ) W l S min > S X,iV+l J W { r B,N+l > r A ,N+l)- 

For the invariant measure of particle iV + 1 we use the the variables £ and W (introduced 
in ref. [|I|) that are defined in the hadron-hadron cm. frame such that the four-momentum 
of particle N + 1 in light cone coordinates (p M = (p + ,p~,~p), p ± = (p° ±p 3 )/V2)) can be 
written as 

In these variables, the soft limit is controlled by £ — > 0, while the limit when particle N + 1 
becomes collinear to particle A is controlled by W— >0 and the invariant measure becomes 

{27Cfx) 2£ d 4 - 2£ p N+1 2 5{p 2 N+1 ) = e~ 2£ d^ (2vr/i) 2£ d 2 - 2£ W (4.43) 

= £W d£diy(sin0)- 2 M0d- 2£ fi, 

where W = |W| and <fi is the azimuthal angle of particle N + 1. In light cone coordinates 
the incoming particles have four-momenta 

Pa = (^yf , 0, o) , p% = (o, a*y~ o) , (4.44) 

hence, £ < xa- The invariant mass squared of particle Q has to be greater then zero, which 
constraints the upper value of £W 2 , ^W 2 < xqXbs/xa with xq = xa — £• For later use we 
record the four- dimensional limit of the integration measure T*™^: 

K^ie = 0) = V N {e = 0)5^ fpg - f>f) , ^d£d^ 2 d0 (4.45) 

V i= i ) | p M =p M +p M„ pM 2 

XW l r min > r min ) W l S min > S X,V+1 J W l r B,V+l > r^JV+lJ- 
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The decomposition of the phase space into initial and final pieces (see eq. ( 4.28 )) naturally 
decomposes the subtraction terms and the finite contribution as well, 



jV 



X[ S oft] m = Yl J [ S ° ft ]m / + E 2>ft]-% 

TTli=A,B 



m f =l 
N 



z[oou] TO = E 2I00UC/+ E 



m,f =1 



A' 



X[fin] = Yl Z[fin] m -f + Yl J [ fin ] 

m/=l rrii=A,B 

In the following subsections, we define these terms precisely. 



(4.46) 
(4.47) 
(4.48) 



4.4 Soft subtractions 

In this subsection we define the X[soft]^ integrals for the cases m = A, B,l, . . . , N and 
x = A, B, 1, . . . , N. We start with the integrals involving the measure X'^-T- 

In the soft limit, the second and third G functions in eq. ( ^4.31| ) take the value one and 
so does the factor (E m + E)/E m in eq. ( |4.32j ). We find 



lim T> 



fin;m , 
N+l 1 



N 



1=1 



V N (e) {2<Kn)- 2s 8 d \p» A +p B -Y:v' 
x(27r/i) 2e d 4 - 2£ p2(5(j9X) 



(4.49) 



m->P 



(N+l;i) r (N+l;f) 



G(r 



x e(min(r ii7V+ i : i = 1, . . . , N, % ^ m) > r m ^ N+l ) 



\mf — >P 



We shall also use the four- dimensional limit of this measure expressed in terms of z and Q 2 : 



hm 2^7(6 = 0) 



(4.50) 



x -dz dQ 2 d^ e(4El(l -z)> Q 2 )9(r 



(N+l;i) ^ JN+l;f), 



> r 

mm ' mm 



x0(min(rj j Ar + i : % — 1, . . . , N, i 7^ m) > r, 



m,N+l, 



This expression is consistent with the soft limit of eq. ( 4.40 ). In the soft limit, z and Q 2 are 
related to the energy and polar angle according to the relations 



z = l-—, Q 2 = 2(1- z)E 2 P (l-cos$). 
E P 



(4.51) 



Using the soft limit of V^^(e), we define the soft terms X[soft]m / (m = A, B, 1, . . . , N) 



as 
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J[soft0 = 

E 

aA, a B, a i>---, a N 



(4.52) 



a. 



N 



V N (e) (2^)- 2£ 6 d b^+^-E^ 



(2vr) 2 

x(27r^) 2£ d 4 - 2e p2 5(^) 



i=l 



ijjf-tF 



6(r 



(N+l;i) (N+VJ) 

' mi n 



xe(min(r iiAr+ i : z = 1, . . . , N, i ^ mf) > r mf)N+1 ) 

x 



S:n,n V®1 Pi PN+l)®( aS mn > s rn,N+l + S n ,N+l] 



y v^ +2 > 

/ j b;mn 
(n=A,B,l,...,N\ 

xL(a A ,aB,XA,x B )S N+ i(p1, . . . 



\mt— >P 



< a< a max < 1/2. 



m y — ».f 
p^=0 



The indices mj and P label the same momentum, therefore, the mj — > P substitution is 
formal in this equation and we shall omit it in the rest of the paper. In eq. ( (4.52p , we 
introduced an upper bound — a times the energy of parton m — for the energy of parton 
iV + 1 in the cm. frame of partons m and n, expressed in terms of invariants. Physical 
quantities will not depend on this bound. The z — 1 — E/Ep > relation introduces an 
upper bound on E which is not present in eq. ( [4.49|) . Thus the use of the soft approximation 
is justified only if 



Q(as mn > s m , N+ i + s n , N+1 )Q(E P > E)S N+1 (d, 



cut) 



(4.53) 



9(aS mn > S mt N+l + S n ,N+l) Sn+1 (^cut). 



for any m and n. When calculating an iV-jet quantity at a fixed value of d CVLt , there always 
exist a finite a max such that this equation is fulfilled if a < a max . The value of a max is easily 
obtained by integrating the measurement function only. 



We now turn to the definition of the soft subtraction terms involving the measure Pjy+i- 
The soft limit, — > is equivalent to taking £ — > 0, in which case eq. ( [4.41| ) takes the form 



A' 



(4.54) 



27T/1 



2; 



d£ dW (sin, 



i=i 



The second line contains the invariant measure of particle iV + 1, so we define the soft terms 
Z[soft]£ (m = A,B,l,...,N) as 
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aA,as,aiv.i]V l Z7 V \ i=l / 

\2e/\4-2e 



X[soft]£ = (4.55) 

/ N 

(2tt)* 

x(27r^) 2£ d 4 - 2£ p2 5(p> M ) 

xe(i^>rrW' w >^ +1 ) 

X E ^S;mn 2) (a)P)P7V+l)0( as mn > Sm,7V+l + S n>JV +l) 

/n=ii,B,l,...,JV"\ 
V n^m ) 

xL(a A ,a B ,x A ,x B )S N+1 (p1, . . . , < a < a max . 

|p/*-»0 

Here, it is useful to choose a max such that 

Q(x A > 0©((^a - Z)x B s > x A ^W 2 ) (4.56) 

xQ(aS mn > S mt N+l + S rii JV'+l)'SjV+l((icut)| 



for any m and n. In this case the generation of the soft phase space in ( 4.69 ) becomes 
simpler. 

The definition of the soft term Xfsoft]^ is analogous. One simply interchanges labels A 
and B in the second G function of the third line of eq. ( |4.55| ). We remark that the sum of the 
N + 2 soft terms Xfsoft]™ 7 and Xfsoft]™ 1 is independent of the phase space decomposition: 



X[soft] m = (4.57) 

/ N 

(2nf 

x(27r^d 4 - 2£ p2 5G»> M ) 

X E ^5;mn 2) (^P,p^ +1 )6(aS mn > S m> jV+l + S n ,AT+l) 



a A ,a B ,a.i r ..,a N J l Z7r J \ i=l / 

\2e j4-2e 



fn=A,B,l,...,N\ 



xL(a A ,a B ,x A , x B )S N+1 . . . , p M ) 

This relation will be the starting point for the evaluation of the integrals over the invariant 
measure of particle iV + 1 in the soft subtraction terms X[soft] m . 

We close the definition of the soft subtraction terms with spelling out the four- dimensional 
limit of the measure (|4.54 ), 

limV™£(e = 0) = V N (e = 0)*« U+p B - f>f) ^dW 2 ^ (4.58) 

x9(r^>rrW^>r^ +1 ) 

This form is consistent with the soft limit of eq. (|4.45| ) . It will be used for giving an explicit 
expression for the finite integral X[fin] in four dimensions. 
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4.5 Collinear subtractions 

In this subsection we define the Xfcoll]^ integrals for the cases m = A, B, 1, . . . , N and 
x — A, B, 1, . . . , N. We start with the integrals involving the measure I^+i and we discuss 
in detail the case X = A. The treatment of case X = B is analogous. 

In the collinear limit p P = zp 1 ^, p M = (1 — z)p^, the first G function in the integration 
measure ( 4.39 ) is zero, therefore, the terms X[coll]^ for x = 1,...,N are defined to be 
zero. Also in this limit, the 6 function Q(rA,N+i > rs,N+i) appearing in the measure V^ 1 ^ 
becomes zero, therefore, the term X[coll];f is defined to be zero. In the same limit the measure 
Djv+i^) becomes 



hmPiSfc) = V N (s)(2n^5 d [zp^+p% - |>f) £ ^ d£ d 2 ^W. (4.59) 



where the momentum fraction z can be expressed in terms of actual integration variables xa 
and £ using the momentum conservation relation for the "+" component of the momenta, 

z = ^ll, (4.60) 

x A 

The collinear pole 2/ sa,n+i equals 2/ (xa£,W 2 ), therefore, the collinear termXfcoll]^ is defined 



as 



I[co\\]i = (4.61) 

/27r/i\ 2e 1 d 2 ~ 2£ W 



E J 7^V N (s)L(aA,aB,XA,x B ) (-?") d £ 



a a, a b, 0.1, ...,ajv 



V £ J x A W 2 



(2nfi)- 2s 5 d (z PA \+p% - f>f) Q(zs5/xa > W 2 )Q(x A > 
E (-ir^+^H;?^ 1 /^ P;A,B,1,...,N + 1) 



X 

a N+l 

\p»=(l-z)p A 
N 



(27r/i)- 2£ 5 rf (pIX + p£ ~ Erf) Q{s5/x A > W 2 )Q(z > 1 



a) 



x ^C{flAl^ 2N) M 

(1 - z) 



|pf=0. 



< 5 < 5 mK * < 1, 



where we introduced a convenient upper bound for the W integral. We have subtracted a 
term in the soft-collinear limit (z — > 1) in order to keep X[coll]^ from having a soft divergence 
when e — > 0. According to eq. Q4.6UD , in this subtracted term the z > 1 — a constraint is 
equivalent to ax a > £ and 1/(1 — z) = x A /C- It is useful to choose 5 max such that 

9({x A - Qxbs > z£s6 max )S N+ i(d cnt ) = S N+1 (d cnt ) , (4.62) 

in which case the generation of the collinear phase space in (f4.69| ) becomes simpler. 
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We shall also use the four-dimensional limit of measure Q4-5SQ 

hmV^i(e = 0) = V N (e = 0)5^ (zfi + f B - f>f ) ^dW 2 d^. (4.63) 

This form is consistent with the collinear limit of eq. ( |4.45 ). It will be used for giving an 
explicit expression for the finite integral X[fin] in four dimensions. 

Next we consider the collinear subtraction terms involving the measure X?™;™ The 
collinear limit of particles m and N + 1 implies taking uJ m ,N+i = uj — > 0, Q 2 — > 0. In this 
limit, in the integration measure (|4.39 ) the first two G functions become one, the third one 
becomes Q(z > 1/2) and the last one becomes S mmf . We find 



lim T> 



fin;m 



;m f 
N+l I 



N 



1/2^ 
2\E P 



V N (e) (2^)- 2 ^ U+^-Ef 
Q 2 x e 



i=l 



5 r , 



(4.64) 



\mt^P 



2r 



[z(l - z)]- £ dz I dQ 2 (sin^)- 2£ dv?d- 2£ fi 

x0(4z(l - z)E 2 P > Q 2 ) G(z > 1/2). 
The collinear pole 2/s mft N+i = 2/Q 2 , so we define the collinear terms T[coll]m / as 
Xfcoll]^ = 



(4.65) 



E 



«i,iiB,ai,...,iiM 



a. 



v 



\ i=i 



(2vr) 2 
xL(a Aj a Bj x Aj x B 



\mt-fP 



/27T/i 



2c- 



2 \ -e 



<3 2 



-2f 



X 



z- £ 0(4z(l - z)£p<5 > Q 2 ) 

x E ^ 2 v + i(^^;A5,i,...,iV + i) 

ajv+i 

x6(z > 1/2)5^+1^, 



^•••,P M ) 



2 ix ^,2x 2 C( 



-6(4(1 - z)E P 8 > Q 



6(01-0)%^,..,/), 



\mr — >P 



(Pm=Pp\ 
\ pM=0 / 



< 5 < <L 



Here we have subtracted the integrand at z — 1 in order to keep X[coll]m / from having a soft 
divergence when e — ■> 0. After making the indicated substitutions, the right hand side does 
not contain the indices m and m/, but the index P. In equation (|4.65| ), for the case 5 = 1 
the upper limit on the Q 2 integration derives from the coscj > — 1 constraint with the use 
of the relations 



Q 2 = 2z(l - z)Ep(l - coscj), 



2(1 - z)E 2 p(l -coscj) 



(4.66) 
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in the collinear and soft-collinear limits respectively. We remind the reader that particle P 
is the massless limit of particle Q. In the collinear limit pt^ = zp P , p M — > (1 — z)pp the first 
G function in the measure ( [4.41|) becomes zero, therefore, the collinear subtraction terms 
X[coll]^ for x = A, B are defined to be zero, therefore, 



jV 



J[coll] m = E 2looU]™/. (4.67) 



mr =1 



Finally we spell out the four- dimensional limit of the measure ( [4.64 ), which is consistent 
with the collinear limit of the measure ( |4.40| ): 



N 

fin;m 



limP^(s = 0) 

1 



T> N {e = 0)6<» U+pS-Erf 



i=l 



(4.68) 



x-dzdQ 2 dcp 0(4z(l - > g 2 )6(^ > 1/2). 
4.6 The finite contribution 

The precise definition of the soft and collinear terms fixes the finite term by eqs. ( |2.3|) and 
( 4.25|) . For the sake of completeness, in this subsection we write it in terms of the actual 
integration variables. The integrand contains at most integrable square-root singularities 
over the integration domain, therefore, it suffices to write the integral in d = 4 dimensions. 
We give the integrals X[fin] mi explicitly for the case = A only. 

X[fin] A = J2 [ 7^2^(0) hdidW 2 d(j)L{a A ,a B ,x A ,x B ) (4.69) 

|e=0 a,A,aB,ai,—,aN l^ 71 *' ^ 



> ) (^+p£-Erf)e(r^>r£. + ^) 

x&(x A > £)e({x A - £)X B 8 > x A £W 2 ) 

x6(£ l) > Sx% +1 )Q(r B:N+1 > r AN+1 ) 

x E * (2N+2) (AB,l,...,N + l)S N+ M,---,P'*N + i) 

a.N+1 

N \ 



r min >r min )^[ r B,N+l > r A,N+l) 

x E ~ 4 r ^ ;c) (a,p)6(as mn > s m , N+ i + s n)N+ i] 

fm,n = A,B,l,...,N\ S m,N+l S n,N+l 

xS N+1 (p1,...,p^) 

| p m=o 

-5 (4) +p£ - £>f) &{z85/x a > W 2 )Q(x A > £) 
y 2 rf iU+/K)+/(«p) 
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x*W3+iO-/ z > P;A,b,i,...,n + i)s N+ M, ..., P n 

U A +jf B - £p?) ®{s5/x A > W 2 )Q{ax A > £) 

s vv \pv=0 ) 

with G [— 7r, 7r] and z given in eq. flOOP . The integral J[fin] s is analogous with the labels 
A and 5 and the momentum components p% + i and p^+i interchanged. The G functions 
assure that the integrand can become singular only when £W 2 — > 0. For the reconstruction 
of the four-momentum see eq. ( |4.42|) . 
The explicit form of the integral X[fin] m ^ is 



lp M =(l-z)f^ 



Z[fin] m ' = E 

|e=0 aA,Hfl,ai,...,ajv 



a 1 

— -^--dz dQ 2 dtp L(a A , a B , x A , x B ) 
{2-K) z 2 



■ 2 \ 2 



/ TV 



.P 5 V^Q, 

x0(s > Q 2 )6(42(l - z)E% > Q 2 



x6(r 



(*f+W > r (N+i;f) 



(N;f) min 
mill ^ ti k,N-\ 



je^ > e 



N+l, 



x0(min(r iiA r + i :i=l,...,N, m f ) > r mf}N+1 ) 
x E * (2N+2) (A,B,1,...,N+ 1)S N+1 

iV 



X 



9(4(1 - z)E% > Q a )e(r& +1 « > r^) 

x6(min(r iiA r + i : z = 1, . . ., N, m f ) > r mftN+1 ) 
2s 



x 



E 



fm,n=A,B,l,...,N\ s m,N+l s n,N+l 

\ m < n / 



-*sr ;c) (ff,p) 



x©(as mn > s mjiV+ i + s n> iv + i)<Siv+i(pi, • • • 



|p^=0 



+0(4z(l - > Q 2 )6(^ > 1/2) 



x- 



Q ajv+i 



iV + 1) 



V = (l-2)p£ 



-6(4(1 - z)Ej,5 > Q 2 )Q{z > 1 - a) 
x2C(a m/ )^— 



(4.70) 



with z G [0, 1] unless otherwise indicated and G [— 7r,7r]. In these equations we used the 
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relation 



E ( 2Smn ' -¥ 2 ^ c \a,p)Q(as mn > s m , N+1 + s n , N+1 ) (4.71) 

( m ,n=A,B,i,...,N\ s m,N+l{ s m,N+l + s n,N+l) 

= E ^mn' C \aiP)®( as mn > Sm,7V+l + S n ,iV+l) 

(m,n=A,B,l,...,N\ S m,N+l s n,N+l 

in order to simplify the sum of the soft subtraction terms. The T>n(0) factor with the 
accompanying momentum conservation S function contains the iV-body phase space, the 
xa-i %b integration (and other trivial factors). This iV-body phase space can be generated 
according to the reader's preference. For instance, one can use the well-known phase space 
generating routine RAMBO J20|. The reconstruction of the momentum of particle (N + 1) 



from z, Q 2 and ip was given in previous subsections (see formulas fl4.35|) , fl4.38|) and (|4.51|) ). 



We have given the precise definition of all terms in eq. ( f4.25|) . The integrals X[soft] and 
X[coll] contain poles in the Laurent expansion in e around zero. According to the factorization 
and Kinoshita-Lee-Nauenberg theorems, these poles cancel against similar poles in the X[2 — > 
N] contribution. In order to see the cancelation of infrared divergences explicitly, we have 
to analyze the integrals over the invariant measure of gluon iV + 1 in the soft and collinear 
contributions, which is the subject of the next two sections. 



5 Soft integrals 



In this section, we evaluate the integrals in X[soft] m (to = A, B, 1, . . . , N) over the invariant 
measure of particle N + 1. 

At the soft point, = the measurement function Sn+i(Pi, ■ ■ ■ ,£> M ) = Sn(Pi, ■ ■ ■ ,Pn)j 
which is the manifestation of the requirement of infrared safe measurement. The only de- 
pendence in eq. ( |4.57| ) on the variables of gluon iV + 1 is in the eikonal factor, therefore, we 
can write X[soft] m as 

2>ft] m = E / l^ V N{e) (27^)-^ (ft +p% - f>f) (5.1) 

a,A,a,B,ai,-,aN \ ~) \ i=l / 

xL(a A ,a B ,x A ,x B )S N (p^,...,p^) E ^™(p)*^ ;c) (a,p), 



(n=A,B,l ] ...,N\ 



where 



Jmn{p)= J(27lfi) 2e d 4 - 2£ PN+1 2 5(p 2 N+1 



(5.2) 



2Sr 



s m,N+l{s m ,N+l + ■Vat+i 



-Q{0iS mn > Sm t N+l + Sn,N+l)- 



This integral is evaluated in the appendix. The result can be obtained exactly, however, for 
our purposes the Laurent expansion in the form, 



JmniP) = 27TC r 



2 

ES, 



— In a In — 

2e 2 e 2e Q 



+ JmniP) + 0(e) 



(5.3) 
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is better suited. The function J mn is independent of e and is very simple: 



2 / 2 °mn 



4 V Ok, 



7i- 

12' 



(5.4) 



Substituting this result for the J mn soft integral into eq. (|5.1|) and using the soft-collinear 
identity, formula ( |4.21]) , we see that X[soft] m assumes very similar form to the X[2 — > iV] 
integral, eq. (|3.10| ): 

X[soft] m = E / ^(^(pSN • • • ,^)(2^)- 2 ^ (p£ + j£ - Erf) (5.5) 

ayi,a s ,ai,...,ajv \ i=l / 



where 



a. 



2 \ e 



xL(a A ,a B ,x A ,x B )—c r . ^ 2 



ES, 



* s m ° ft (a,p), 



^ ft (a,p) = * {27V) (a,p) 



+ E « ;c) («,p) 



V n^m } 



-2C(a m ) In a 



1 j / ^mn 



(5.6) 



6 Collinear integrals 

In this section, we evaluate the integrals in X[coll]™ (m = A, B, 1, ... , N) over the invari- 
ant measure of parton N + 1. Before going into the details, we have to make a remark. 
The collinear subtraction terms were defined using the four-dimensional expressions for the 
collinear limit of the squared matrix element. That was sufficient for the evaluation of the 
X[fin] integral. Strictly speaking however, the subtraction scheme applied in this paper is de- 



fined in d dimensions. It was shown in ref. ]T5| that with making use of process independent 
transition terms, one can use four- dimensional expressions for the helicity independent part 
for the collinear limit of the squared matrix element except for the Altarelli-Parisi splitting 
functions, P a /b that have to be calculated in d dimensions. As for the helicity dependent 
part, the analysis of its general structure in d dimensions shows it has vanishing azimuthal 
integral in d dimensions ||, Therefore, we drop the helicity dependent part of ^c^mn^ 
in the following considerations. This causes some inconsistency in our notation, but the 
physical cross section remains unchanged. 

We start with the evaluation of the integrals in X[coll]™ (m = A, B) over the invariant 
measure of parton N + 1. At the collinear point W— > with £ fixed, the measurement 



function Sn+i(Pi- 



We change integration variables in the first 



term of eq. ( |4.61| ) form xa (which is hidden in T>n and in the luminosity factor) and £ to 
xp = xa — £, and z with z defined in eq. ( |4.60|) . The Jacobian for this transformation is 
xa/z. The lower limit for the z integral is defined by the xa < 1 relation, hence z > xp, 
while the upper limit is obviously one. The limits on xp are zero and one, just as was on 
xa- After this change of variables we can rename the index P to A (and simultaneously 
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the flavor index a A to b). In the second term of eq. Q4.61 ), we change variable from £ to z. 
Keeping the helicity independent term, we can now write X[coll]^ as 

(6.1) 



Z[ooU]i 



a A ,a B ,ai,...,a N J Z7r \ i=\ J 



1 dz 




[( z XA ) 







-2e 



x(27r/i) 2e / tj79 9 (s5z 2 /x A > W' 



W 2 



- ((1 - z)x A ) 2e Q(z > 1 - a)L(a A , a B , x A , x B ) 



r H 2— 2e W 

x(2n^j '±-^L®(s8/x A >W'< 



2C(a A ) 
1 - z 



Evaluation of the integral over W (see eq. ( [A. 12 )) results in 



Z[ooU]5 



(6.2) 



E / V N (e)SM, . . . y N )(2n^5 d U+p B - Erf) 

!,ai,...,Ojv \ i=l / 

P aA /b(z, e)L (b, a B , x B ^j 
2C(a A ) 



l0(,>x,)E^ 
z 2 ^ u(a A ) 



-<d(z > 1 — a)L(a A , a B , x A , x B )- 



1 - z 



In order that we could combine this contribution with the collinear factorization counter 
term for hadron A, we use the relations 



(47T) 



er(l-e) 



cr 



2 

ES, 



and 



X £ f 1 dz(l-z) 
Jo 



-2e 



i + .n (Sk) + O(e) 



2C 



(6.3) 



/(*)P(z, e)e(« > z) - /(!)— 6(z >l-a) 



(6.4) 



d*/(z) 



1- z [l-z]+ 



-/(l)[ 7 + 2Clna] 



+e< / dz 



+2C/(1) [lnX(ln(l -x) - In a) + ln 2 a-ln 2 (l - x)} \ + 0{e 
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The latter equation holds if (1 — z)P(z, 0) — > 2C as z — > 1 and for any function /(z) that is 
not singular in z — 1. The function is defined by the relation P(z,e) = P(z) +eP'(z). 

For the complete collinear integral X[coll]^ we obtain 



J[coll]^ 



(6.5) 



E / V N (e)S N W, ■ ■ ■ MQw)-*? U+p B - Erf 



x \ 7T- C r 



2 \ £ 
^ \ ,t,co11/^ 



2tt V.Q es 



2 , (a,p)L(a A ,a B ,x A ,x B ) 

ES/ 



2 \ e 



2tt^ \Q 2 



^ {2N) (a, p)L(a A , a B , x A , x B ) 



X 



f 1 dz\nX A (z)(±j: 

Jx A \Z b 



u(b) L(b,a B ,x A /z,x B ) p _ 2C(a A ) \ 

u(a A ) L(a A ,a B ,x A ,x B ) aA b 1-z J 



Jxa z 2 ^ u(a A ) L(a A ,a B ,x A ,x B ) aA/b 
-2C(a A ) (\nX A (0)\n^— + ln 2 {l-x A ) -In 2 a) | + 0(e), 

where X A (z) = /U 2 /((l — z^rr^s^) and 

n n &P) =^ (2W Ha,p)^[7M+2C(a A )lna]+ln^^ [ 7 (a A ) + 2C(a A ) In a]) . (6.6) 

When two final state particles become collinear, pi^ = zp P , p^ = (1 — z)p P , the mea- 
surement function S^+iipi, ■ ■ -P 41 ) = Sn(Pi, ■ ■ ■ ,Pn), ■ Keeping the helicity independent 



term in the collinear integral, we have 
Z[coll]™ = 

E / \v N (e)S N (pt, . . .,p%)(2nn)- 2e 5 d (p A +p B - Erf 

V ap,...,aj\r / 



(6.7) 



i=i 



\m^P 



x 



x 



(27T) 



-L(a A ,a B ,x A ,x B ) 



m->P 



Q 2 



X 



0(4^(1 - z)E 2 P 5 > Q 2 )Q(z > l/2)J2P am/ap (z,e) 



-6(4(1 - z)E 2 P 5 > Q 2 ) Q(z > 1 - a) 
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The Q 2 , if integrations can be calculated immediately. The required integral is given in the 
appendix, eq. ( |A.9| ). The remaining integral over z is also straightforward. After performing 
these steps and leaving out the formal m — >• P substitutions, we see that the term X[coll]™ 
has a form very similar to that of X[2 — > N]: 



J[coU]; 



(6.8) 



aA,a,B,ai,...,a,N \ i=l 



xL(a A , a B ,x A ,x B )—c T 



ft 



2 \ £ 



2tt^ \Ql s , 



where 



ES 



e \4E15. 



Z am (a)¥ 2N \a,p) + 0(e) 



(6.9) 



with Z am (a) given in eq. ( |A.10| ). 

At this point we see that the sum of the integrals X[2 — > N], X[soft] and X[coll] is free of 
any poles of e, therefore, it can be calculated in d = 4 dimensions: 



(X[2 -> AT] + X[soft] +X[coll]) 



(6.10) 



e=0 



E / 2M°)S*(P?, ■ ■ ■ '^)^ (4) (jtf +P& - Erf) L(aA,a B ,x A ,^)vI/( 2Ar )(a,p) 



X 1 + 



a, / * 



2ti\ tf( 2JV )(a,p) 



fm,n=A,B,l,...,.fV\ 



^( 27V )(a,p) 



+ E [7(«m) + 2C(a m ) In a 



m=l 



hi 



+ E [7(0 + 2C(a m ) In a 



AE15 

ln^fe 
/i 2 



iv r 



E 



m=l . 



7T 



2C(a m ) In 2 a + — + 7 '(a, 



+ E 2C K 

- f 1 d^lnX A (^) 

Jx A 

- f 1 dz\nX B (z) 



lnX m (0)ln 



1—Xr 



+ In (1 — x m ) — In a 



1 x . L(b,a B ,x A /z,x B ) ~ 2C{a B ) 

2 2 ^ L(a A , a B ,x A ,x B ) 1 ~ z 

1 x - L(a B ,b,x A ,x B /z) ~ 2C(a B ) 

z z L(a A , a B ,x A ,x B ) 1-z 



1 A 1 sr o>b,xa/z,xb) p, ( s 
x A 2 b ijya A , a B , x A , x B ) 

_ r 1 , 1 \ - L(a B ,b,XA,x B /z) 5, , " 
l s % 2 V L(a A ,a s ,a; A ,x B ) aB ' b{Z} / 
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The logarithmic dependence on the unphysical parameters a and 5 in this equation gets 
canceled when this contribution is combined with the integral Z[fm] (eq. ( [4.48 )) in order to 
obtain the infrared safe physical cross section at next-to-leading order : 



a = J[fin] + ( J[2 -> iV] + J[soft] + J[coll] ) 



(6.11) 



e=0 



e=0 



We can also demonstrate the independence of the iV-body integral of the auxiliary pa- 
rameter Q|s explicitly by making use of eqs. ( |3.12| ), ( [4.21 ), ( |5.4| ) as well as the definition of 
the £2 function, formula (13.91): 



(J[2 -»■ N] + J[soft] + J[coll]) 



(6.12) 



e=0 



E / V N (0)S N W,...,p%)6^ U+p B -j: P ?)L(a A ,a B ,x A ,x B )¥ 2N \a,£) 

aA> a B,ai,---,a>N \ i=l / 



X 1 + 



o:. s 



2vr 



,27V 



/V 



f) 

4tt/ 



+ E E 

m=l ( n=A,B,l,...,. 

+ E E 



hel col 



+ In a In 



m=A,B (n=A,B,l,...,N\ L 



^-Q(smn) + In a In 



/V 



+ E 



lM^—- 5 + T -C(a m ) — 



*( 2A 0(a,p) 

Y(Om) - 7(Om) 



E 

m=A,B 



2C(a m ) (ln(l-x m ) In Xm ^ ^ m " >sS _ ZL^ _ 7 ( flr 
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- / dzlnX A (z) 



X A 
1 

1 

XA 



dz hxX B (z) 



E 

6 

■E 



L(b,a B ,x A /z,x B ) 5 



L{aA,a B ,XA,XB) 1 



2C(a A ) 



L(a B ,b, x A ,x B /z) ~ 
L(a A ,a B ,x A ,x B ) 



a B /b 



Z - 



Z 

2C(a B ) 



^E 



L(b,a B ,x A /z,x B ) ~, 



L(aA,a B ,XA,x B ) 



K A M 



i dz l_ L { a B ,b,x A ,x B ,z) p \\ 
x B z z L(a A ,aB,XA,x B ) Bl )) 

This equation together with eqs. ( |4.69|) and (f4.70|) define explicitly those integrals that are 
needed for the calculation of a jet cross section at next-to-leading order in perturbative QCD. 



7 Numerical results 



In this section we present some numerical results of the first non-trivial application of our 
algorithm, namely the calculation of three-jet cross sections in e + e~ annihilation. Thus our 
results can be compared with those of ref. [ f2~l|j . 
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Figure 1: Coefficient of (a s /(27r)) 2 for the thrust, C-parameter distributions and the y cnt 
distributions of the Jade E and fc_i_-clustermg algorithms. The thrust distribution is multi- 
plied by (1 — t), the C-parameter distribution is multiplied by C, while the distributions of 
the clustering algorithms are multiplied by y cnt . The dotted histograms show the statistical 



errors. 



We use the matrix elements of ref. [2^] for the construction of the various \I/ functions. 
The algorithm can easily be altered for performing jet cross section calculations in the case 
of e + e~ annihilation. One simply drops all terms in the integrals (|4.70|) and (|6.12| ) that 
carry A or B indices, leaves out the xa, %b integrations from V^(0) and then the sum of 
integrals ( 4.70 ) and ( 6.12Q immediately gives the physical cross section. We implemented 
such an algorithm in a Monte Carlo program. The results are in good agreement with those 
of ref. As an example we show the next-to-leading order coefficients for the thrust, 



C-parameter distributions multiplied by (1 — t) and C respectively and distributions for the 
Jade E and k± jet clustering algorithms multiplied by the jet resolution parameter in fig. 1. 
We find that the numerical convergence is similar to the program of ref. 
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8 Conclusion 



In this paper we have presented a simple generalization of the subtraction method of ref. 
H for the calculation of any infrared safe physical cross section in perturbative QCD. The 
apparent conflict between the need for important sampling in order to achieve sufficient 
numerical precision and the increasing difficulty of performing partial fractioning in the 
tree-level next-to-leading order matrix element was overcome by a decomposition of the 
phase space such that in one region only one Lorentz invariant of the external momenta 
can become singular in the calculation of an iV-jet observable. We wish to emphasize the 
simplicity of our algorithm: the necessary analytic integrals are rather trivial, while the 
numerical implementation is only a little more complicated than a tree-level Monte Carlo 
program. 

We have given all the necessary integrals that define any next-to-leading order QCD jet 
cross section explicitly (see eqs. ( |4.69| ), (|4.70|) and (|6.12|) ). The phase space integrations 



in these integrals can be programmed for any number of jets. Once having such a master 
program the only ingredients that have to be changed in a modular fashion are the 

• the Born-level and next-to-leading order tree matrix elements in four dimensions (^I^ 2 ™) 
(n = N, N + 1) functions); 

• the color linked Born matrix elements (\|/( 2Ar;c ) functions); 

• the non-singular part of the one- loop helicity amplitudes (.4ns functions); 

• the S n (n = N,N + 1) measurement functions. 

The algorithm can be trivially altered for calculating QCD jet cross sections in other 
processes, like in e + e~ annihilation or deep-inelastic scattering. As an example, we have 
shown results of the next-to-leading order thrust, C-parameter and jet distributions in the 
case of e + e~ annihilation. The numerical convergence was found to be similar to that of the 



program of ref. in the case of e + e~ annihilation into three jets. In the case of three-jet 
production in hadron collisions such a benchmark calculation does not exist yet. In order to 
demonstrate the applicability of the algorithm in such calculations, in a companion paper 



231 , we give results of a next-to-leading order calculation of three-jet cross section in hadron 
collisions in the simplified case of pure Yang-Mills theory. The structure of the algorithm is 
essentially the same when quarks are included therefore, the conclusions are expected to be 
similar in the full QCD case. 
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Appendix 

This appendix is added for the readers' convenience. It contains a collection of the 
integrals that were used in the main text. The explicit evaluation of these integrals is rather 
trivial, therefore, minimal details are given. 

First we calculate the soft integral 

J mn (p) = J (27t f i) 2£ d 4 ~ 2£ p N+1 2 5(p 2 N+1 ) (A.l) 

2s mn . . 

X -, ; rk3{aS mn > S m ,AT+l + S n>N+ i). 

S m ,N+l{S m ,N+l + S n ,N+l) 

for m, n — 1, . . . , N, m ^ n. The integrand depends only on Lorentz invariants therefore, 
the integral can be calculated in any frame. We choose the "m-n" system, where the four- 
momenta of gluon m, n and N +1 take the form (first four components are energy, z, x and 
y components of the three-momentum) 



^ = ^(1,1,0,0,...), (A.2) 



# = ^(1,-1,0,0,...), (A.3) 

p% + i = E(l, cosfl, sin d cos </?, sin^sin^, . . .), (A. 4) 
where the dots mean d — 4 zeros. With this choice the invariants are 

S m ,N+l = y/s^E(l - cos$), (A. 5) 

y/s^E(l + coati). (A.6) 
Consequently, in the "m-n" system the integral takes the following simple form: 

J mn = J (^) 2£ ^ (sin^)- 2£ d^ (siny)- 2 ^ d^ i ^ tf e(aV^/2 > E). (A.7) 
This integral is easily evaluated and one obtains the exact result: 

Jm " - ? U^J * rtl^i)' (A ' 8) 

In the main text, we have also used the following integrals: 

/ (wY > ^ = -T Q ™W=7)> (A - 9) 



\E 2 P ) Q 2 v ^ max ^ ' v ^ r £ ^ max r(i -£)' 

= / ^ (i - ^ 2£ £ A/a(*)e(* > i/2) - ^^e(^ > 1 



-7(a) - 2C(a) In a + £7' (a) + 2£C(a) In 2 a + y J , 



— a) 

z 

2\ 



(A.10) 
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where 



and 



The derivation of these results is sufficiently simple so that we can omit the details. 
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